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Abstract. We discuss a simple analytical model of the turbulent boundary layer 
(TBL) over flat plane. The model offers an analytical description of the profiles of 
mean velocity and turbulent activity in the entire boundary region, from the viscous 
sub-layer, through the buffer layer further into the log-law turbulent region. In 
contrast to various existing interpolation formulas the model allows one to generalize 
the description of simple TBL of a Newtonian fluid for more complicated flows of 
turbulent suspensions laden with heavy particles, bubbles, long-chain polymers, to 
' include the gravity acceleration, etc. 
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Various problems of environmental and engineering hydrodynamics 
call for a simple analytical model for the TBL over flat plane, that can 
adequately describe from a unified viewpoint the mean velocity and 
turbulent activity profiles in the entire boundary layer. In this paper 
we analyze in details such a model, announced in [1] in connection 
with a problem of drag reduction in dilute polymeric solutions. The 
model is based on the balance equations for mechanical momentum 
and kinetic energy. Aiming maximum possible simplicity of the model 
. we neglect the spacial energy transfer in favor of the energy production 

and dissipation. This makes our model local from a physical viewpoint 
and algebraic from an analytical side. We also suggested a closure 
that links the Reynolds stress with the density of kinetic energy. In 
contrast to known interpolations between the resulting formulas for 
the mean velocity profile in the viscous and turbulent sublayers we 
suggest a uniform model for the rate of energy dissipation at the point 
of the formulation of the model. Besides the physical transparency, our 
approach allows straightforward generalization of the model for more 
complicated flows of turbulent suspensions laden with heavy particles, 
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bubbles, or long-chain polymers, inclusion of the gravity acceleration, 
etc. 

The basic Navier-Stokes equation (NSE) for the fluid velocity U(r, t) 
can be written as 



= -Vp + fiV 2 U, (1) 



where p is the fluid (air) density, p = p(r, t) - the pressure and p is the 
dynamical viscosity. In this paper we follow the standard strategy of 
Reynolds, considering velocity as a sum of its average (over time) and 
a fluctuating part: 

U(r, t) = V(r) + u(r, t) , V(z) = (U(r, t)) . 

The objects that enter our model in the planar geometry are the mean 
shear S(z), the Reynolds stress W(z) and the kinetic energy K(z); 
these are defined respectively as 

S(z) = , W(z) = -p(u x u z ) , K(z) = P -{\u?) . (2) 

Here x, y, and z are (horizontal) streamwise, spanwise, and (vertical) 
wall-normal directions. 

Integrating the stationary NSE for the mean velocity V(z), one gets 
a well known exact relation [2], that describes the point-wise balance 
of the flux of mechanical momentum: 



fiS(z) + W(z) = V{z) . (3) 

In the RHS of this equation we see the total flux of the mechanical 
momentum V; in the LHS we have the Reynolds stress and the viscous 
contribution to the momentum flux. At large Reynolds numbers one 
can usually neglect near the surface the production of Viz) due to the 
pressure gradient or by some other reasons. If so, 

V(z) =V = V(0) . (4) 

Having in mind a lower part of the planetary boundary layer (PBL), 
we consider Eq. (4) as a boundary condition at the ground level z = 
instead of a given value of a free stream velocity at the upper boundary 
of PBL. The value of Vq gives natural "wall units" u T , r and l T for the 
velocity, the time and the length: 
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Introducing so-called "wall normalized" dimensionless objects 

,+ s f, v = ^, s+ s ^, w^^L, k^JL, 

(5) 

we can rewrite Eq. (3) as 

S + {z + ) + W + {z+) = 1 . (6) 

A second relation between S(z), W(z) and K(z) is obtained from 
the "point-wise" energy balance: 

e prod = ^dis , locality approximation. (7) 

in which we neglected the spacial energy transfer term, e^. The detailed 
analysis, see, e.g. Fig. 3 in Ref. [3], shows that in the log-law turbulent 
region this term is small with respect to the energy dissipation term 
e dis : ^tr ^5 O.ledis- Clearly, in the viscous sub-layer the mean velocity 
is fully determined by the viscous term and thus the influence of the 
energy transfer term can be neglected. For simplicity of the model we 
will neglect e tr term also in the buffer layer (where the ratio e proc i/edis 
is between 1 and 1.8). We will show below that the locality approxima- 
tion (7) does not essentially affect the resulting mean velocity profile 
and Reynolds stress. 

The energy production rate e pro d in Eq. (7) describes the energy 
flux from the mean shear flow to the turbulent subsystem. In the plane 
geometry it has a simple (and well known) form that follows from the 
NSE (1): 

eprod = W(z)S(z) . (8) 

It is also well known that the kinetic energy K dissipates due to 
viscosity at the rate 

edis(^) = 

In the viscous sub-layer the velocity field is rather smooth, the gradient 
exists and thus can be reasonably estimated via the distance to the 
surface as 1/z. In other words, in this region we can write 

Cdis (9) 

where 

K{z), v=- p , (10) 

with a being some dimensionless phenomenological constant of the 
order of unity. 
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In the buffer sublayer and in log-law turbulent region the energy 
cascades down scales and is finally dissipated at the Kolmogorov (inner) 
scale that is much smaller than the distance z. Therefore the contribu- 
tion to the energy dissipations from all scales, smaller than z, is equal 
to the energy flux, which we denote as efl ux . This means that outside 
the viscous sublayer Eq. (9) has to be supplemented by an additional 
term, e flux : 

e dis (z) =eZ(z) +e flux (z) . (11) 

Notice, that in the buffer sublayer, both contributions in (11) to £&\ s {z) 
are equally important, while in the log-law turbulent region the direct 
dissipation of energy of turbulent eddies of the largest scale z in the 
system, given by Eq. (10), is negligibly small with respect to the non- 
linear energy flux efl ux . Clearly, in the viscous sublayer, Eq. (11) also 
should work, because the nonlinear contribution, ea ux (z) is negligibly 
small with respect to the linear one, e^,(z). We believe that Eq. (11) 
is more than just interpolation for the energy dissipation between the 
viscous sublayer and log-law turbulent region. As we will show below, 
the model (11) gives an uniformly reasonable description of the rate of 
the energy dissipation in the entire boundary layer. 

To make this description constructive one has to evaluate in Eq. (11) 
the energy flux efl ux (z). This can be done by standard Kolmogorov 
1941-type dimensional reasoning: 



^^m^jMJm. (12) 

t(z) z V P 

Here r{z) is the typical eddy turnover time at the height z equal to 
the turnover time of the largest eddies (of scale z) at this height and b 
is another dimensionless constant of the order of unity. Thus Eqs. (8, 
10, 11) and (12) allows us to rewrite the energy balance Eq. (7) in the 
entire boundary layer as follows: 



W(z)S(z) 



K(z) 



K{z) . (13) 



P 

In the dimensionless "wall-normalized" objects (5) this equation reads: 



W + (z + )S + (z + ) 



K + (z + ) . (14) 



Now we have two balance Eqs. (6) and (14) for three objects, S + (z + ), 
W + (z + ) and K + (z + ). Two of them, W + (z + ) and K + (z + ), are different 
components of the same Reynolds stress tensor {ufu^). Therefore it is 
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naturally to expect that in the scale invariant region (which in our prob- 
lem is the log-law turbulent region) these objects will have the same z 
dependence and thus their ratio will be z-independent (dimensionless) 
constant: 

- <V) =► cl ■ (15) 

Notice that this ratio is bounded from above, c 2 (z + ) < 1, by the 
Cauchy-Schwarz inequality. In fact the expectation (15) with ~ 0.28 
in the log-law turbulent region is in a good agreement with numerous 
laboratory and nature experiments, see e.g. book [2] and with many 
DNS data, see for instance below Fig. 3, taken from Ref. [3]. 

Needless to say, that various Reynolds-stress based closure proce- 
dures lead to the same result, Cqo = const, in the log-law turbulent 
region, in which Cqo is expressed via yet another phenomenological 
constants. In our simple model we prefer to use Eq. (15) as a basic 
closure. Moreover, we argue that we can safely use Eq. (15) not only 
in the log-law turbulent region, where it is definitely valid, but also 
in the buffer layer and even in the viscous sublayer, where Eq. (15) 
is violated. The reason is simple: the larger the deviation of the ratio 
W + (z + )/K + (z + ) from the constant c^, the less important become 
relation (15) itself in the momentum and energy balances. Below in 
this paper we account numerically for the real z + -dependence of the 
ratio W + (z + ) / K + (z + ) and demonstrate that this is insignificant for 
the mean velocity and the Reynolds stress profiles. 

Equation (15) allows us to represent our model, Eqs. (6) and (14), 
in terms of just two unknowns, S + and K + or W + . We choose W + 
instead of K + , because the Reynolds stress is responsible for the tur- 
bulent transport of the mechanical momentum and thus plays a more 
important role in the wall bounded turbulence than the kinetic energy. 
Notice, that in the wall turbulence W + is positive definite, since the 
momentum flux is directed toward the surface. 

In terms of S + and W + Eqs. (6) and (14) read: 

1 = S + (z+) + W + (z + ), (16) 







clS + (z+)-r + (z + ) W + (z + ) . (17) 



Here T + (z + ) can be considered as an effective damping rate of the 
turbulent fluctuations: 

r+(z+) = (^) 2 + -^-^+(3+), (18) 

and c 2 XJ S + (z + ) clearly represent the energy influx rate. 
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The basic equations of our model (16) and (17) have two solutions: 
a laminar and a turbulent one. In the laminar solution there are no 
turbulent fluctuations: 

W + (z + ) = K+(z + ) = 0, (19) 
S + (z + ) = 1, V + (z + ) = z + : laminar solution. 

The stability condition with respect to appearance of the turbulent 
fluctuations requires that the damping rate, r + (z + ), at a zero level of 
turbulence is larger (or equal) than the pumping rate, c 2 ^ [recall, that 
in the laminar solution S + (z + ) = 1]: 

(~+^J — c oo > stability of the laminar solution. (20) 

This equation shows that the laminar solution (19) is stable near the 
surface, for z + < z+, where in our model 

Zy S = — is the upper boundary of the viscous-sublayer. (21) 

Coo 

Recall, that the energy transfer is neglected in the model. Therefore 
it is not surprising that Eq. (19) demonstrate no turbulent activity in 
this sublayer. 

There is however a turbulent activity in the rest of the boundary 
layer: 

W + (z + ) > , for z + > z^ s , mixing layer & log-law region, (22) 
in which Eq. (17) gives 

cls + (z + ) = r + (z + ). 

This relation together with definition (18) yield: 

b 



2 

C 2 



Dividing this equation by c 2 ^, and using definition (21) for z+, one 
finally gets instead of (16, 17) a new set of coupled equations: 

1 = S + (z + ) + W + (z + ), (23) 
^(* + )= f4V + ^p. (24) 

\ Z + KZ + 



Here we introduced another dimensionless parameter k 

k=&, (25) 
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which in our model is nothing but the von-Karman constant, that 
defines the slope of the logarithmic mean velocity profile in the log-law 
turbulent region. Notice, that the final system of coupled equations (23) 
and (24) have a minimum possible number (just two) of phenomeno- 
logical constants, z+ and k (that are some combinations of initially 
introduced three parameters, a, b and Cqo). Indeed, any models of wall 
bounded turbulence have at least two phenomenological parameters, 
see,e.g. [2]. For example, the famous "logarithmic law of the wall" 



V + {z + ) = k- 1 In z + + B, for z 
k « 0.436, Ba 6.13, 



>30, 



(26) 



contains the von-Karman constant k and the intercept B, with exper- 
imental values in Eq. (26) taken from [4]. 

Let us show, that unlike (26), Eqs. (23, 24) describe the velocity 
profile in the entire boundary layer and not only in the log-law turbu- 
lent region. Eliminating W + from Eqs. (23, 24) one gets a quadratic 
equation for S + with two solutions. The physical one has the form: 



2k 2 (z+) 2 - 1 + J4«2 



+ 1 



(27) 



2k 2 z+ 2 

Now Eq. (23) immediately gives an expression for W + , which is valid 

for z + > z+ : 



2k 2 



W + 



(4 S Y 



+ i 



Uk 2 



y.+- 



-OvVf 



+i 



2k 2 z~ 



(28) 



One sees that at z + = 
the laminar solution (19): S + (z+ s ) = 1, VF + (; 



the turbulent solution (27, 28) coincides with 



0, as expected. To 



get the mean velocity profile, we integrate Eq. (27) matching the result 



with the laminar solution (19) at f 



. Fortunately, the expression 



for the mean shear (27) allows analytical integration. In is convenient 
to present the result of this integration in the form, similar to the 
logarithmic law of the wall (26): 

V + (z + ) = k- 1 In Z(z + )+B- A(z+), for z + > z+ . (29) 

Here functions Z(z + ), A(z + ) and intercept B are given by 

1 



Z(z+) 



A(z+) 



2 

2 k 2 z: 



+ 

-2 



- & + {^r 

+ Ak[Z(z+) -z+]-l 



2k 2 z+ 

B = 2z+ - k- 1 In [e (1 + 2 kz+)/& 



for z 



for z 



oo , 
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Note that Eq. (29) pertain to the whole z + > z+ domain, meaning 
both mixing sublayer and log-law turbulent region. By taking the ex- 
perimental values of k and B we compute z^ s ~ 6 to be compared with 
the experimental value of 5.5 db 0.5, cf. [2]. The resulting mean velocity 
profile for the entire boundary layer, Eqs. (19) and (29), is shown in 
Fig. 1 as the solid line. The excellent agreement with the experimental 
and numerical data in the entire region of z + indicates that our balance 
equations are sufficiently accurate. 

In Fig. 2 we show analytical profiles of the Reynolds stress W + (z + ), 
Eq. (28), in comparison with DNS. Due to the limited value of the fric- 
tion Reynolds number in the DNS data [3], ReA = 590, we normalized 
the Reynolds stress using the local value of the momentum flux 

^ + '="°( 1 -£t)- 
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The same normalization was used in Fig. 5 for the kinetic energy. Notice 
that the type of normalization [with Vq or V(z + )] does not affect the 
mean profile (due to its slow, logarithmic dependence on z + ). There- 
fore, in the DNS data for mean velocity in Fig. 1 we used the simple 
normalization with Vq. One can see in Fig. 2 an excellent agreement of 
our analytical results with the DNS results. The minor discrepancy is 
observed only in the viscous sublayer, z + < z+ s : in our simple approach 
the Reynolds stress and the turbulent kinetic energy are identically 
zero in this region, see Eq. (19). This stems from the disregard of the 
energy transfer term in the energy balance equation (7), which gives 
a non-zero level of turbulent activity close to the surface. As one sees 
however from Figs. 1 and 2, the account of the spacial transfer term 
does not lead to a considerable value of the Reynolds stress in the 
viscous sublayer and does not affect the mean velocity at all. 

Another assumption that fails in the viscous and buffer layers is 
the approximation of constancy of the correlation coefficient = 
W + /K + . Note, however, that the expressions for the mean shear (27) 
and the Reynolds stress (28) remain valid even for z + -dependent cor- 
relation coefficient Cqo c(z + ). In this case z+ and k should be 
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understood as z + -dependent functions: 

_j- _|_ / & / _|_\ c (z~^~) 

In fact, we have chosen c(z + ) = const only to make possible an analytic 
expression for the mean velocity profile, Eq. (29). In our model we 
can easily account for the "realistic" z + -dependence of the ratio c 2 
integrating Eq. (27) numerically. 

The actual dependence c(z + ), shown in Fig. 3 by blue dashed line, is 
taken from the public available statistical database, produced in Ref. [3] 
by DNS of the NSE for high-Reynolds turbulent flow in the channel 
geometry. One sees that c(z + ) decreases toward the surface. This fact 
can be understood by a series expansion for VF + (z + ) and K + (z + ) for 
z + —>■ (see, e.g. [2]). This expansion shows that near the surface 
W + (z + ) and K + (z + ) behave as 

W + (z + ) ~ (z + f, K + {z + )~{z + f, 

and therefore 

c 2 (z + ) ~ z + , near the surface. (30) 

An origin of these dependencies is quite simple. At the surface, 
the rms values of the horizontal projections of the turbulent veloc- 
ity \J (uf 2 ) and \J (uf 2 ) are zero according to the no-slip boundary 
conditions and grow with the height as z + . This is not the case for the 



vertical projection, y (uf ). The incompressibility constraint dictates 
that 

du z ( du x du v s 



dz V dx dy 
and therefore the vertical projection increases only as z 2 : 



yA^~^f>~* + , ^f)~(z+f. 

Accordingly, 

W+{z +) = _ (u+ u +) „ ^ {u + 2 ) (uf 2 ) ~ (z+f , 

while the kinetic energy is dominated by the horizontal turbulent ve- 
locities, 

K + (z+)^ l -[(uf) + (uf)] ~(z+) 2 , 
in agreement with Eq. (30). 
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Figure 3. Blue dashed line - DNS data [3] for c 2 (z + ), red dash-dotted line - sug- 
gested fit Eq. (31) for c 2 (z + ), and horizontal black solid line - asymptotical value 
=0.28 [3]. 



To demonstrate how the actual dependence c 2 (z + ) influence the 
mean velocity and turbulent activity profiles we suggest the following 
interpolation formula, 



1 — exp 



(31) 



that has just two parameters, asymptotic value c^, for z + — ► oo and the 
slope c^/z^ of the linear dependence (30) for z + — » 0, near the surface. 
Dependence (31) with = 0.28 and z~^. = 24 is shown in Fig. 3 by red 
dash-dotted line. One sees that Eq. (31) closely fits the DNS data in the 
entire boundary layer and thus can be used as a realistic representation 
of the W + /K + ratio in our model to find the analytical representation 
for the improved profiles of the Reynolds stress and kinetic energy, 
Eqs. (15, 28) and then, after numerical integration of Eq. (27) to get 
the improved mean velocity profile. 

The comparison of the resulting profiles is given in Fig. 2 (Reynolds 
stress), Fig. 4 (mean velocity) and Fig 5 (kinetic energy). The profiles in 
the simple model with c(z + ) => Cqo are denoted by black solid lines, the 
improved profiles with z + -dependent coefficient c(z + ) - by red dash- 
dotted lines and the DNS profiles - by blue dashed lines. One sees in 
Fig. 4 that all mean velocity profiles nearly collapse: our approximation 
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has no effect on the function V + (z + ). For this object we prefer to take 
c(z + ) = Coo and to have a fully analytic model. The profile of the 
Reynolds stress, as one sees in Fig. 2, is improved in the viscous sublayer 
and is affected very little in the rest of boundary layer by account for 
z + -dependence of c(z + ). As for the profile of the kinetic energy, Fig. 5, 
the neglect of the actual z + -dependence of the coefficient c(z + ) leads to 
a significant underestimate of the kinetic energy in the buffer sublayer. 
In particular, the simple model does not exhibit a peak of K + (z + ) in 
this region. If this peak is essential for some particular reasons, one 
should account for the actual z + -dependence of c{z + ) at the expense 
of simplicity. 

We have to stress, that the effect of the spacial energy transfer, 
neglected in our simple approach, is absolutely insignificant for the 
mean velocity profile (see Figs. 1 and 4). It has a minor importance for 
the profile of the Reynolds stress (Fig. 2) in the buffer layer and plays 
a more important role for the profile of kinetic energy, decreasing the 
amplitude of its peak and increasing the value of K + (z + ) from the left 
of its maximum, in the viscous sublayer. Therefore, the necessity to 
account for the transfer term should be evaluated for each particular 
problem in hand. 
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Figure 5. Mean profiles of the kinetic energy K + (z + ): black solid line - analytical 
profile with c = c x red dash-dotted line accounts for the z + -dependence (31) of c 
and blue line - is the DNS profile [3]. 

In conclusion: In this paper we discussed in details a simple model 
of a turbulent flow over a flat plane that offer an analytical description 
of the mean velocity, Reynolds stress and kinetic energy profiles in 
the entire boundary layer. The calculated profiles exhibit an excel- 
lent agreement with the results of the laboratory experiments [4] and 
DNS [3] of the NSE. We discussed the effect of the approximations 
made in the analytical description of the profiles. We found a simple 
functional form for the experimental z + -dependence of the correla- 
tion coefficient c(z + ), that allows to relax the approximation of the 
constancy of this coefficient. The profiles calculated with the help of 
this z + -dependent coefficient c{z + ) account for all physically important 
features of all three profiles in the entire boundary layer region. The 
physical transparency and simplicity of the model allow its generaliza- 
tion for turbulently flowing suspensions, as it was demonstrated on the 
example of the problem of drag reduction by polymers in Refs. [1, 5, 6]. 
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